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deriving results appropriate to the strict limit \a\ — > oo in coordinate space. We then extend these 
results to finite a using once-subtracted momentum-space integral equations. We also discuss the 
implications of our results for experiments that probe three-body recombination in Bose-Einstein 
condensates near a Feshbach resonance. 
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I. INTRODUCTION 



Systems in which the two-body scattering length, a, is much larger than the range, I, of 
the underlying interaction are said to be close to the "unitary limit", \a\ — > oo. In that limit 
the two-body scattering saturates the bound set by the requirement that the S-matrix be 
unitary. The dynamics of such systems for energies E ~ h 2 /{Ma 2 ) (M is the particle mass 
is universal, in the sense that it is independent of any details of the underlying two-body 
interaction. One obvious example of this "universality" is that for a > all such systems 
have a two-body bound state (dimer) with binding energy E D « h 2 /{Ma 2 ). 

A non-relativistic system of three identical particles with \a\ 3> I displays a rich set of 
universal features pQ. In the limit \a\ — > oo the system has an infinite tower of three-body 
bound states (trimers). Their binding energies are: 

B n = (e- 2 ™T~ n *£n* (1) 

where B n * denotes the energy of an arbitrary state in the spectrum, and the numerical 
constant So ~ 1.00624. The nature of the geometric spectrum Q is a consequence of the 
discrete scale invariance of the three-body problem in the joint limit \a\ — > oo, I — > 0. Efimov, 
who made this discovery, pointed out that this discrete scale invariance is also relevant for 
finite a and that it should be observable in systems for which l/\a\ ^ 0, but is still <C 1 
Hi- 

Evidence for Efimov physics in cold atoms was recently presented by the Innsbruck group 
[I] . Using a magnetic field to control the scattering length, they measured the recombination 
rate of cold 133 Cs atoms in the lowest hyperfine state and found resonant enhancement at 
a ~ — 850ao, which can be attributed to the existence of an Efimov trimer at the three- 
body threshold. Approaches which employ a zero-range model and therefore account for 
the separation of scales between a and I by taking the limit / — ► are very successful in 
describing this experimental data [5J E] . But a systematic calculation of the corrections to 
these results due to finite-range effects is desirable, in order to identify the limitations of 
these previous calculations and also to extend the range of applicability of approaches based 
on the scale hierarchy I <C \a\. 

In recent years Efimov's results have been rederived in the framework of an effective field 
theory (EFT) [7]. This EFT employs the ratio l/\a\ as a small expansion parameter and 
allows for a systematically improvable and model-independent calculation of observables 
of few-body systems with a large scattering length. The leading order (LO) of this EFT 
reproduces all universal features previously discovered by Efimov. 

In this paper we will exploit effective theory methods to derive the correction linear in 
the effective range to this universal behavior ([I]). In the EFT expansion this corresponds to 
the next-to-leading order (NLO) piece. The value of the two-body effective range, r s , sets 
the size of these corrections. We evaluate their impact perturbatively, and so obtain for the 
three-body bound-state spectrum a result of the form: 



B n = B 



II* 



— J + r s G n ( — ) + 0[(«* r s 
ft* / V ft* 



(2) 



where «;* = \/MB nJf is the binding momentum of the n*th state in the three-body bound- 
state spectrum, 7 = 1/a (at LO) is the binding momentum of the two-body bound state 
that is present for a > 0, and the functions F n obey pQ 

F„(e s o x ) = e s o F n _ m (x). (3) 
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We will show that the functions G n obey approximate relations similar to Eq. ([3]), since 
discrete scale invariance connects the G n 's for different values of n. Knowledge of the 
functions F nJf and G n therefore provides information on the behavior of the bound-state 
spectrum in the vicinity of a narrow Feshbach resonance to a precision given by (r s /a) 2 . In 
this work we will give results for the function G n * from the unitary limit up to the value of 
a where the three-body bound state under consideration hits the threshold. In particular, 
we will show that this function is zero at particular values of a. 

We do this by formally studying systems in which the effective range, r s , is significantly 
larger than the range I of the underlying two-body potential, i.e. the scale hierarchy is 
/ <C \r s \ <C \a\. It is important to point out that in such a three-scale problem it is possible 
to be in the unitary limit, but still have range effects present in the result, i.e. to not have 
achieved what is sometimes called the "scaling limit", where r s — > and I — > 0. In fact, 
in the EFT formulation of the problem the range of the underlying two-body potential, I, 
provides a limiting distance, because the EFT does not specify any details of the two-body 
potential, and so it does not accurately describe the three-body dynamics in those regions of 
configuration space where any two-body distance is less than /. Hence it is natural to place 
a cutoff A on the momentum-space integrals present in the EFT calculations, and choose 
A ~ 1/7. In our case we are interested in systems with I <C \r 8 \, and so we will ultimately 
take the limit A — > oo = I — > in order to obtain our results. The effective range, r s , will, 
however, be kept finite. 

The inclusion and effects of range corrections in the EFT with contact interactions have 
been studied extensively. In [S] range corrections in the three-nucleon system were calculated 
perturbatively in momentum space. A formalism which employs a partial resummation of 
range effects was developed in Refs. [HI QUI EH] and observables were calculated up to next- 
to-next-to-leading order (N2LO) in the l/a expansion. Hyperradial coordinates and the 
Wilsonian renormalization group were used to re-derive the leading-order amplitude and 
develop a power counting for sub-leading effects in Ref . [12] . 

In Sec. [IT] we will briefly review the hyperradial formalism and show analytically that in 
the limit of infinite scattering length no linear range correction exists, i.e. that G n (0) = 
for all n. In Sec. Ill we will review a formalism that allows us to compute numerical results 
for this linear correction away from the unitary limit and present numerical results for the 
three-body bound-state spectrum and the function G n . In Sec. IV we will lay out how the 
formalism developed for the analysis of Sec. Ill can be extended to the computation of range 
corrections for the three-body recombination rate in a cold gas of identical bosons. We will 
end with a summary and discussion of our findings. 



II. NEAR A FESHBACH RESONANCE: RESULTS IN THE UNITARY LIMIT, 
WITH NON-ZERO EFFECTIVE RANGE 

In this section we first review the use of hyperradial coordinates for the three-body 



problem in the unitary limit (Subsection II A). Our presentation closely follows that of 



Ref. PQ. In order to fix notation and methodology for the subsequent subsections we derive 
the potential between the single particle and the pair in hyperspherical coordinates, and show 
that it behaves like one over the hyperradius, R, squared for |r s | <C R <C \a\. We explain 
why — in the joint limits r s — > and \a\ — > oo — this results in the Efimov spectrum described 
above, and compute the wave functions in this limit. A somewhat different explanation of 
these results — one based on Efimov's original work [21 13] — can be found in Appendix A of 
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Ref. [12]. 

We then show how to compute perturbative changes to this spectrum as we move away 
from the scaling limit. In Subsection II B we derive the terms in the hyperspherical potential 
that scale as r s /R 3 , and compute the matrix elements of this potential between the wave 
functions obtained with the pure 1/R 2 potential. These matrix elements diverge in the 
ultraviolet, but we show that if the spectrum is renormalized by demanding that short- 
distance physics at R ~ r s remains unchanged when we move from LO to NLO then there 
is no change in the Efimov spectrum at 0(r s ). 



A. Hyperradial formalism: a review 

The material of this subsection is based on the discussion of the hyperradial formalism 
in the recent review of Ref. p] . For three particles of equal mass the Jacobi coordinates are 
defined as: 

T id = Yi- Yj ■ r k>ij = r k - -(n + r,-) , (4) 

where the triple (ijk) is a cyclic permutation of the particle indices (123). The hyperradius 
R and hyperangle a k are then defined by 

R 2 = \{* 2 i2 + 4 + 4) = \v% + ; a k = arctan [ ^4 ) . (5) 



3 via ■ ■ -ai, 2 -ij -r 3 Mi' ~» I 2 \r k 



1.1 1 



In the center-of-mass system the Schrodinger equation in hyperspherical coordinates is given 
by 

Tr + T ak + + V(R, Si)) *(R, a, SI) = E*(R, a, Si), (6) 



2MR 2 



with 



2M V dR 2 AR 2 

h 2 1 / d 2 \ 
Ta = 2MR 2 ^[-d^- A ) Sin2 ^ (8) 

A M / = + -4^. ( 9 ) 

sin a k cos z a k 

where SI = (0ij,(f>ij,8k,ij,<j>k,ij) an d the Ls that appear in Eq. Q are the usual angular- 
momentum operators with respect to these angles. 

Assuming that the potential V depends only on the magnitude of the inter-particle sep- 
aration we write 

V(r h r 2 , r 3 ) = V(r 12 ) + V{r 23 ) + V(r 31 ). (10) 

We now employ the usual Faddeev decomposition of ip for three identical bosons and neglect 
subsystem angular momentum 

V(R, a, SI) = ^(R, a x ) + ^{R, a 2 ) + ^(R, a 3 ). (11) 

Putting (|ll]) and Q into @ we get 

(T R + T ai - E)^(R,a 1 )+V(V2Rsina 1 )(ip(R,a 1 )+^(R,a 2 )+^(R,a 3 )) =0. (12) 
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Here we have chosen, for definiteness, k = 1. 

We can now exploit the fact that ip{R, a{) is independent of f 2 3 and f 12 3 to obtain the 
simplified Faddeev equation 



(T R + T a — E)ip(R, a) = -V{V2R sin a) [ip(R,a) + 



V3 J\i 



. 7T I 7T 



!f- Q| sin 2a' 



sin 2a 



•0(i?, a') da' 



(13) 

The solution of this equation can be expanded in a set of eigenfunctions of the hyperangular 
operator, i.e. 



ip(R, a) 

where the functions 4> n satisfy 
d 2 



1 



R 5 / 2 sin(2a 



J2fn(R)MR, 



(14) 



da' 



- Xn(R) 



2M R 2 

(R,a) = --^-V(V2R sin a) 

1% 



71" I 7T „ I 

4 /-2-|6- Q l 

l (i2,a) + -^= / <p n (R,a')da' 

V3 J|f-a| 

(15)" 

with boundary conditions <p n (R,0) = <f> n (R, f) = 0. Meanwhile the hyperradial functions 
/„(i?) satisfy 



a 2 is 



+ 



+ V n {R) 



fn(R) + £ 



_2M V <9# 2 AR 2 
with the hyperradial potential V n (R) defined by 



d 

2P nm (i?)— + Q nm {R) 



f m (R)=Ef n {R), 
(16) 



K(i?) = (A„(i?)-4) 



2MR 2 ' 



(17) 



and P nm (R) and Q n m(R) potentials that induce coupling between different hyperradial chan- 
nels p]: 



Pnm \R) 
Qnm (R) 



2M 



J2 G nl(R) f' 2 ''da<fil(R,a)^Ll> m (R,a 

u J 

2 



n 2 



L, JO 



2M 



with 



Grim \R) 



n/2 



da (j)* n (R, a)4> m (R, a). 



(18) 
(19) 

(20) 



Equations (18) and (19) make it clear that P nm and Q nm vanish exactly if the n 's are 



independent of R. As we shall now show, if the potential V is short ranged, then the n 's 
are indeed independent of R 1 . 



1 Below we show that P nm (R) and Q n m(R) contribute to the energies of the Efimov spectrum at order 

{Ts/Rf- 
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For hyperradii R which are much larger than the range over which V is non-zero, say Z, 
we can consider a large enough that V(y/2R sin a) is zero. The solution of Eq. (|15j) in this 
region is then 

a) « sin yXn[ a) . (21) 



ligh) 



Note that we have implicitly assumed that A n becomes independent of R, an assumption 
which will be justified a postieri in the limit R <C \a\. 

On the other hand, as long as we choose R^$> I, there is always a region a < l/{y/2R) <C 1. 
Here h 2 \ n (R) / (2M R 2 ) is much less than V(y/2R sin a), which typically will be of order 
h 2 /(2M1 2 ). Consequently for this a domain the differential equation for <j) n {pt) becomes 



h 2 d 2 tf ow) 

2MR 2 da 2 



+ V(V2R, 



a 



ow) 



--^V(V2Ra)<> 



Oil 



(high) 



7T 



(22) 



Here we have used the fact that a <C 1 in this region to approximate the integral in Eq. ( 15 ), 
and to make the replacement sin a — > a. If the derivatives are all of natural size, then the 
corrections to this equation are suppressed by a 2 . 



The homogenous version of Eq. (|22|) has the solution 

A^ k=Q {V2Ra) 



low) . 



a 



(23) 



where ipk(f) 1S the wave function of the two-body system for two-body energy and Ax 
is a constant. Adding a particular solution of the inhomogeneous equation we find 



ow) 



a) 



Aip (V2Ra) ^sin (\f\n 

v 3 ^ 



(24) 



We now consider Eq. (22) for values of a such that l/R < a <C 1. In this domain the analysis 

re now beyond the rang 

cos(kr) + cot S sin (At)] . 



that led to Eq. (24) still applies. But because we are now beyond the range of V, we have 

sin(Ar + 5(h)) sm5(k) 



A; 



A: 



(25) 



As k — > this yields ipo{ r ) = r — a, and we can use this asymptotic two-body wave function 
in Eq. (24). This gives 



ow) 



a) 



A(V2R 



so 



a — a) — sin 



7T 



(26) 



But, since V = in this region, this result must be consistent with Eq. (21). This is 
achieved by the choice 



.4 



sin 



7T 



2\ ' 



(27) 



which ensures that 0„(a) is continuous across the boundary between "low" and "high" 
solutions at a ~ l/R, and the condition 



cos 



v/3a; 



: sin 



A "6 




A, 



sin 



x 4) -• 

2/ a 



(28) 



on A n , which ensures that <f) n (a) has a continuous first derivative as a — > l/R. We note 
that if these equations are satisfied A n , and hence <fr n , is independent of R for R <C \a\, as 
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promised. Indeed, as long as Eqs. (28) and (27) are satisfied the form (21) is the result for 
(p for all a such that a > l/R. Solving Eq. (28) in the limit R <C \a\ we find the lowest 
eigenvalue 



Ac 



si[l + 1.897 : 



(29) 



with so = 1.00624.... This is the only negative eigenvalue, and hence only this channel 
potential is attractive. So, if we now focus on the unitary limit, where \a\ — > oo, we have 
A = 



Sq. The hyperradial equation (16) in this channel then becomes 



d 2 



2M V OR 2 



4 + 

R 2 



f (R) = Ef (R). 



(30) 



This equation will hold for R I. If we desire a solution for negative E the necessity to 
have /o be normalizable mandates that: 



/o (0) (*) 



RK, 



IS 



2kR), 



(31) 



where the superscript (0) indicates that we are working in the unitary limit, while the 
subscript refers to the solution for the hyperchannel corresponding to Ao, which is the only 
one that supports bound states. The binding energy of these bound states is related to the 



K of Eq. (31) by 



2,^2 



B 



(32) 



Since the attractive l/R 2 potential produces a spectrum that is unbounded from below 
some other short-distance physics is needed in order to stabilize the system. If the two-body 
potential is known this short-distance physics is provided by the two-body potential V, that 
becomes operative for R ~ I. But an alternative approach is to add an additional term to 



Eq. (30) that summarizes the impact of the two-body V. Here we take this potential to be 



a surface delta function at a radius 1 /A 



V SR (R) = H (A)A 2 5 R 



(33) 



with Hq adjusted as a function of A such that the binding energy of a particular state, 
say B n * (with a corresponding given by (32)), is reproduced. Note that since Vsr is 
operative only at distances R ~ 1/A it corresponds to a three-body force. (See Ref. [7] for 
a realization of this in a momentum-space formalism.) 

In physical terms we anticipate I ~ 1/A, since we know that once we consider hyperradii 
of order 1/A the potential V starts to affect the solutions. Our goal here is to derive uni- 
versal results, that are independent of details of V. Consequently we will use the procedure 
described in the previous paragraph to fix the value of H for a given A, and then examine 
the residual A-dependence of our results. If we can take the limit A — ► oo without it signif- 
icantly impacting our calculations of observables (formally, if all effects scale with negative 
powers of A) then these observables are not sensitive to the physics that takes place in the 
three-body system for R ~ I. Although we know that our theory is not valid in the region 
R ~ I, its predictions that have only small (~ 1/A) corrections due to that short-distance 
physics will be "universal" and not significantly affected by the inaccuracies in the short- 
distance region. As we shall see, the excited-state spectrum is one such observable — which 
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is not surprising, since these states have very small binding energy with respect to the scale 
I. 

Given that our focus is on predictions of the theory that are independent of details of 
V we can consider the extreme case and take the limit Z — > 0. In this limit the form of 
K iso as R — > guarantees that once H is fixed to give a bound state at \E 3 \ = B w , the 
other binding energies in this hyperradial eigenchannel form a geometric spectrum. Namely, 
B n = h 2 K 2 n /M with 

K n =(e-^ S0 ) n ~ nt ^, (34) 

with n* the index of the bound state corresponding to «*. Thus Eq. is obtained. Eq. (34) 
will hold for all K n such that K n <C A. (Note that now the subscript on k denotes the index 
of the bound state in adiabatic channel zero.) The continuous scale invariance of the 1/R 2 
potential has been broken down to a discrete scale invariance by the imposition of particular 
short-distance physics on the problem through the short- distance potential (33) |13j . 



B. The Linear Range Correction in the Unitary Limit 

The solutions derived in the previous section are valid in the strict limit \a\ — *■ oo. In 
this section we consider the corrections to this limit that are order (r s /c*). In performing 
this computation we shall assume that the range of the two-body potential I is zero 2 . The 
calculation we are performing therefore corresponds to the limit where I <C |r s | <C \a\. Since 
the natural scale for the cutoff A ~ l/l we will consider A|r s | 3> 1. 

The first-order correction to EFT results in the \a\ — > oo limit was already discussed in 
Ref. [8]. However, in that work Hammer and Mehen considered r s /\a\ effects, and so the 
first-order corrections they discuss go to zero in the limit \a\ — > oo that we are considering 
here. As we shall show here, in principle there is a first-order correction proportional to 
r s K* that could survive if one considers \a\ — > oo but keeps r s finite. However, it turns out 
that the discrete scale invariance that relates wave functions for different bound states at 
r s = guarantees that the particular operator in question has zero matrix element in the 
limit A — > oo. 



To prove this we first observe that in the limit A|r s | ^> 1 the form (21) applies for all 



values of a, as long as Eq. (28) is obeyed. There are then two possible sources of corrections 



to the Efimov spectrum (34) that are linear in r s . First, we must reconsider Eq. (22), and 
take account of the previously neglected term h 2 X n /(2MR 2 ). As we shall see this leads to 
the conclusion that 4> n (a) is not independent of R, and so then we must also check that P nm 
and Q nm remain zero at the level of accuracy we work to here. 



Reconsidering Eq. (15) we assume that V ~ h /{2ml ) and consider the region l/R 

2 V(V2Rsma)^<P^ (r, 



a C 1. We can then write 
d 2 . 2MR 2 



- A n H —V(V2Rsma) 



da 2 ' h 



n 2 v 'vr n v '3. 

(35) 



2 



One might be concerned that this precludes consideration of a positive effective range, because of the 
bound derived in Refs. [TJJ[TS]. But that bound applies only to energy-independent potentials. If the 
underlying two-body potential V(R) is energy dependent positive effective ranges can be generated no 
matter how small I is. 
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The approximate form adopted here for the integral on the right-hand side of Eq. (15) is 
itself accurate to relative order a 2 , meaning that the right-hand side omits only terms of 
overall order a 3 . Therefore the first correction to Eq. (35) is 0(l 3 /R 3 ) relative to leading. 
This equation thus goes one order in a ~ l/R beyond Eq. (22). 

In general the homogeneous form of Eq. (35) is difficult to solve. But, for the specific 
case of I — > 0, we need only consider the effects of V for a — > 0. There sin a = a still holds. 
Hence the solution of the homogenous form of Eq. (35) for any a > becomes: 

^ ow \R, a) = A{R)i/; kn (V2Ra), (36) 
/x 



with ipkir) defined by Eq. (25) and k n = ^=j|. Meanwhile, the particular solution of Eq. (35) 
is unchanged up to the accuracy to which we work here, and so 



(37) 



4™\R,a) = A(R)^ kn (V2Ra) - sin ( y/X^j + 0(a 3 ). 



V3 ~ V v "-"6. 

Substituting in the form (25) for ipk n ( r ) we see that 0n° w ' 1 will still match smoothly with the 
large-a solution (21), provided that 

A{R) 



sin (5 (k n )) 



sin 



An: 



7T 



and 



cos 



7T 



v/3a; 



: sin 



7T 



-R k n cot(5(k n )) 
R A„.r 




4 R 



+ 



A 
R 3 



sm 



7T 




A, 



An" 



sm 



7T 



(38) 

(39) 
.(40) 



In Eq. (40) we have assumed only that the coefficients of all terms in the effective-range 
expansion after the 0(k 2 ) one scale with I. Under this assumption the scattering length and 
effective range determine the eigenvalue A„, up to corrections suppressed by (l/R) 3 . 

In this case Eq. ( |40| ) can be solved by noting that X n (R) has a self-consistent solution in 
the region \r s \ <C R <C \a\ 



An( J R) = Af(0) ^l + 7n ^+^ n 



R 



+ ... 



(41) 



with Xn\R) is the solution of Eq. (28), which here we need only at R = 0. By a perturbative 
evaluation of Eq. (40) in powers of 



r, we find that: 



Ai 0) (0) 7 r, 
4 



encodes the 0(r s ) correction to A n in the limit \a 
eigenvalue is 

2 



A„( R) = -4 ( I + 1.897^ + 0.480^ + 



(42) 

oo. Hence, the value of the lowest 

(43) 
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where we have omitted corrections of second order in R/a and r s /R and have not considered 
r s /a corrections. Similarly, the second-lowest eigenvalue can be determined as: 



Ai(.R) = 19.94 1 



0.1872- + 0.9333^ + -- 
a R 



(44) 



In what follows, we will focus only on the lowest eigenvalue Xq(R). 
With Aq(-R) in hand we can compute 



[a;R) 



sm 



a 



(45) 



for a > l/\a\, up to corrections of order (r s /R) , with Ao given by Eq. (43). For later 
purposes it makes sense to normalize 0o(a; R) such that [16J 



(a; R) da 



irrespective of the value of R. Since the overall normalization of 
done, resulting in: 



(46) 



is arbitrary this is easily 



4> (a;R) 



N(R) 



sm 



7T 



— a 



with 



N 2 (R) 



71 



JO ITT 



(47) 



(48) 



where jo is the spherical Bessel function of zeroth order. 

The fact that A n now depends on R opens the possibility that the hyperchannel coupling 
potentials P nm and Q nm will no longer be zero. Clearly both can be at most O (^), and so 
at leading order in r s we need only consider P nn (R) and Q nn (R)- Given that our focus here 
is on the bound-state spectrum this means we only have to examine Poo(R) and Qqo(R). The 
use of normalized <f)o(a;R) can then easily be seen to guarantee that P 00 (R) is identically 
zero [16J. A brief computation also shows that 



Q 0Q (R) =0 + O 



R 2 



(49) 



and so we need not concern ourselves further with these adiabatic coupling potentials in 
pursuing our calculation of effects linear in r s . 

Consequently, in the unitary (\a\ — > oo), but not scaling (r s ^ 0) limit the hyperradial 



Eq.(16) becomes 



1L 

2M 



d 2 s 2 + l 4£ r, 



OR 2 



R 2 



R 3 



f (1 \R) = E {1) f m (R) 



(50) 



where the superscript one indicates that we have worked only to first order in r s . The form 
of Eq. (50), but not the specific pre-factor of the 1/R 3 potential, was proposed by Efimov 



on dimensional grounds |17j . 

The momentum-space evaluation of first-order (in r s /\a\) effects in Ref. [8] neglected the 
operator corresponding to the 1/R 3 potential that appears in Eq. (50). The authors of 



Ref. [8] argued that the piece of the two-body scattering amplitude that survives in the 
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limit \a\ — > oo and is proportional to r s does not have a two-body bound-state pole and so 
cannot lead to effects in the three-body system. As we shall now see the conclusion that the 
1/R 3 potential does not affect the spectrum in first-order perturbation theory in the unitary 
limit is correct, although it is not justified to neglect this piece of the three-body dynamics 
entirely, as was done in Ref. [8]. 

Now — to the extent that the r s /R 3 term is a perturbation — the shift in the Efimov spec- 
trum near a Feshbach resonance (i.e. near |o| = oo) can be evaluated by computing the 
matrix elements of this perturbing 1/R 3 potential between leading-order hyperradial func- 



tions (31). In other words, we have, for the first-order shift in the bound-state energy: 



2M 



ti 2 



A5« = s 2 r s £ j 



dRti 0) \R)^, 



(51) 



with EP = -(^ + ABk 1 ') the energy of the nth bound state up to first order in r s . The 
f!? ] that appears in Eq. (plh is the normalized zeroth-order hyperradial wave function of 



the corresponding bound state [TJ: 



fW(R) = 2U S ^^R^K lS0 (V2n n R). 



TCS 



(52) 



Unfortunately the integral in Eq. (51) is divergent. This divergence can, however, be 
absorbed by a modification of the short-distance potential (33) to include terms ~ r s 



we now add to Vsr(R) of Eq. (33) a piece: 



V^(R)=H 1 (A)A 2 5[R 



A 



I.e. 



(53) 



with Hi ~ r s . 



We include this potential in Eq. (50), and then treat it in perturbation theory. Per- 



turbation theory is only manifestly valid if |r a |A <C 1. But, as we will show below, after 
renormalization we can take the limit A — > oo in our perturbation-theory expressions. Al- 
though the 1/-R 3 potential and Vg^ cannot individually be treated in perturbation theory 
in this limit their combined effects remain perturbative for arbitrarily large A. Indeed, we 
will show that the A dependence in the perturbation-theory calculation is confined to terms 
of 0(1/ A) and so it is always small. The complete first-order-perturbation-theory result is 
therefore independent of details of short-distance dynamics, i.e. it is "universal" in the sense 
we are using the term here. 

At finite A the total first-order shift in the energy of the nth bound state in the zeroth 
eigenchannel is then: 



2M 



A 



^l A 2 f (0) 2 

2 Jn V A 



(54) 



where the theory has been regulated at a distance 1/A and we shall seek to remove this cutoff 
(i.e. take A — > oo) at the end of the calculation. We shall, however, calculate throughout in 
the limit that A ^> (indeed, that A is much larger than the binding momentum of any 
states we are interested in). The dimensionless hi in Eq. (54) is defined as: 



hi 



AHiM 

h 2 sl r s £ 



(55) 
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and the subscripts on the hyperradial wave functions indicate which binding energy B n they 
correspond to. From now on we will drop the superscript on these / n 's, since for the rest of 
this Subsection all wave functions are to be evaluated at leading order. 



Substituting for /„ from Eq. (52), and using the short-distance form of K iSQ (x): 

x 



K iso {x) 



71 



with [I] 

we get: 
2MAB { n 1] 

h 2 



a 



S0 sinh (7 r S0 ) Si H S0l H7lJ +a0)+O( ^ 



11 f T(l + is ) 

-s ln2 - -arg — — 

2 2 \r(l — is Q ] 



(56) 



(57) 



s sinh(7rs ) 



7T 



dR- 



KfjV2K n R) h 1 (A)A . 



R 2 



} sin [s \n (-^ ) + n„ 



(5« 

The result of the integral inside the square brackets can be expressed in terms of Hyper- 
geometric functions. Since we are only interested in the result for A ^> n n we expand these 
functions in powers of K n /A and keep the first two terms to obtain: 

r°° dr 1 
I{z) = / ™Kl Q (x) = -U(z) + zh{z) + 0(z% (59) 

with the functions J_i and I\ encoding non-analytic dependence of the coefficients on z: 



h(z) 



7T 



2s sinh(7rs ) 



7T 







^=^= cos ( 2sq hi(z/ a/2) + 2a + arctan.(2so)") 
l + 4s£ V / 



(60) 



4s (l + Sq) sinh(7rs ) 



x cos I 2s ln(z/v2) + 2a — arctan 



3s 



l-2ag 



(61) 



Inserting Eq. ( 59 ) into the expression for the energy shift we have: 
2MAB { n 1] 



h 1 



s sinh(7rs ) ( ^ ( V2 



7T 



A 



lKl h (V2Kn 



A 



A 



/ii(A)A 



o sin [sq In ( — ) + a 



(62) 



It follows from Eqs. (62) and (60) that the coefficient of the linear divergence can be 



forced to zero if we choose: 

MA) 1 



sin (s lnx + « ) 



1 - 



a/1 + 4^ 



cos (2sq ln(a;o) + 2a + arctan(2so)) 



(63) 



12 



with xq = k^/K. With this choice of hi we can write: 



2MAB { n } 

h 2 



(64) 



where the function g is proportional to I\(\/2K n / &) up to small corrections and so is gener- 
ically of order one — although it does have log-periodic dependence on K n /A. Similarly the 
function / inherits the log periodicity of 



f(ze°°) = f{z), 
while the renormalization condition guarantees: 

f(x ) = 0. 



(n* — n)7r 

Since, at leading order, Kyi — s ° , Eqs. (|65|) and (j66h result in: 



0. 



(65) 



(66) 



(67) 



Therefore once the linearly divergent part of the shift is renormalized away for one state, 
it will — at first order in perturbation theory — be zero for all bound states in the Efimov 
spectrum. Consequently, 



4£o?V< 
A 



-g 



X 



(68) 



Numerical results generated for the first-order shift to the spectrum in the unitary limit 



conform to this pattern [18]. Eq. (68) holds as long as we renormalize in such a way as to 
keep fixed. 

Therefore a renormalization procedure that leaves one of the Efimov states unperturbed 
by r s corrections leaves all states unperturbed in the limit A — > oo. This conclusion is a 
consequence of the discrete scale invariance of the Efimov spectrum. Because the leading- 
order spectrum has discrete scale invariance, and states are related to one another by a scale 
transformation of a factor e~ n n *^*o, the integrals (51) for different states n are also related 
to one another by that factor, i.e. 



e *o ABl 1] . 



(69) 



Thus, if ABn] is forced — by our renormalization procedure — to equal zero, we must also 
have 

AB^ = for all n. (70) 

In order for this conclusion to hold it is critical to use a regularization and renormalization 
prescription which respects discrete scale invariance. The cutoff-plus-delta-function proce- 
dure we adopted does this in the limit A — > oo. At any finite A the violations of discrete 



scale invariance result in 1/A-suppressed corrections to Eq. (70). However, even if A is kept 
~ 1/7, the perturbative corrections to the spectrum would only be ~ rj^, and so is of 
higher order than we are considering here. 
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III. THE LINEAR RANGE CORRECTION FOR ARBITRARY SCATTERING 
LENGTH 



In the previous section we focused on the limit in which \a\ —>■ oo, i.e. ft* ^> l/\a\. 
The hyperradial formalism is particularly well-suited to the analysis of that limit, as the 
hyperradial potential is a power law (or, more generally, sum of powers of r s /R) for all values 
of R. In this section we use a momentum-space formalism [20] to obtain the corrections to 
the Efimov spectrum for arbitrary values of 1/ (an*). Numerical calculation of scattering and 
bound-state observables is straightforward, and it is easy to compute range corrections — 



regardless of the value of a [ID] . We review the relevant formalism in Section III A Results 



for the range correction are presented in Sec. HI B and analyzed for their content in terms 



of universality in Sec. Ill C 



A technical point is that to do this we employ the effective-range expansion around the 
bound state pole 

A;coa = -7 + 7^ 2 + 7 2 ) + --- , (71) 



rather than the one around k = 0, since Eq. (71) facilitates the treatment of the two-body 



threshold cut generated by the three-body scattering equations. The parameter p is equal to 
the effective range r s up to corrections ~ r 2 /a, which means that setting p = r s is adequate 
for our NLO calculation. Therefore in what follows we do not distinguish between p and r s . 



A. Subtracted Momentum Space Equations: A Review 

For zero-range, I = 0, two-body potentials the usual momentum-space Faddeev equations 
for the atom-dimer scattering amplitude reduce to the equation first derived by Skorniakov 
and Ter-Martirosian (STM) [19] 3 : 

A^(q,q';E) = Z(q,q';E) + t dq" q" 2 Z(q,q" , E) ^^(l ^°\q", q';E) , 

(72) 

where A^ denotes the leading-order atom-dimer scattering amplitude, Ed is the dimer 
binding energy, A a momentum cutoff which regularizes the integral equation, e a positive 
infinitesimal, and 

Z^^-^J ^/^-fX (73) 
qq'n \q + q — qq — I 

Finally, in Eq. (72), S m (E; q") = S (0) (e - denotes the function: 

with Ed = h 2 j 2 /M and E = h 2 k 2 /M. Note that here we have Ed (equivalently \a\) finite, 
in contrast to the analysis of the previous section. Note also that if 7 < then the pole 



The STM form appears different, but is equivalent [TIED], 
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corresponding to the dimer is on the second sheet of the complex energy plane, and so does 
not correspond to an actual bound state, but instead to a "quasi-bound state". Regardless 



of the sign of 7 Eq. ( 72 ) may be made more compact by introducing the dimer propagator 



K J E + E D + ie K ' 

Although now 50 years old the STM equation has recently been re-derived as the leading- 
order result for the atom-dimer scattering amplitude in an effective field theory (EFT) with 
contact interactions alone. This EFT is the appropriate low-energy theory for non-relativistic 
systems with a ^> I and k( < 1, where k denotes the momentum scale of observable of 
interest [7J. 



The EFT analysis of Ref. [7] showed that Eq. (72) leads to strongly cutoff-dependent 
results for the bound-state spectrum. In coordinate space this is reflected by the necessity 
to fix the phase of the hyperradial wave function at short distances before a prediction about 
the three-body bound-state spectrum can be made. As in the analysis of Sec. [XXJ it transpires 
that fixing a single three-body observable (in what follows we choose the threshold atom- 
dimer scattering amplitude) is sufficient to remove the cutoff dependence of the spectrum. 
We will now show how to achieve this using a once-subtracted version of the integral equation 



(72). Our presentation follows Ref. [20] . 
We begin by considering: 

am r A 

A<%, 0; -E D ) = Z(q, 0; -E D )-— J dq" Z(q, q"; -E D )S^(-E D , ; q")A^\q" , 0; -E D ) , 

(76) 

We can insert the required three-body input by demanding that the scattering amplitude 
reproduces the correct particle-dimer scattering length at threshold. The on-shell amplitude 



at E — —Ed obeys Eq. (76) with q — 0, but the value of A^(q, 0; — Ed) is fixed if the 



atom-dimer scattering length 03 is taken as input. Therefore by subtracting the equation at 



q = from Eq. ( 76 ) and relating to a 3 we obtain the subtracted equation at threshold: 



A<®{qA-E D ) = 3 ^ + A[Z](q,0;-E D ) 

o'y 



3 

Here 



am r A 

' ' dq"A[Z](q } q";-E D )S^(-E D ;q")A ( - \q",0;-E D ).(77) 







A[Z](q, q'; E) = Z{q, q'; E) - Z(0, q'; E), (78) 

and denotes the atom-dimer scattering length. 

Now we can determine the full-off-shell amplitude at threshold by exploiting the symme- 
tries of A. Before subtraction the full off-shell amplitude at threshold satisfies 

A<%, q'; -E D ) = Z(q, q f ; -E D )-— J dq" Z(q, q"; -E D )S^{-E D] q")A^\q" , q'; -E D ). 

(79) 

Since Z(q,q';E) = Z(q',q;E), we also have „4 (0) (0,g; -E D ) = A (0) {q,0; -E D ). Using this 



fact as input, together with the solution of Eq. (77), produces the subtracted version of 
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Eq. ((79]): 

A W (q,q';-E D ) = A^(0, q'; -E D ) + A[Z](q, q'; -E D ) 



4M rA 



o 



3 

which determines the fully-off-shell amplitude at threshold 



dq"A[Z](q, q"; -E D )S®(-E D ; q")A^(q", q'; -E D ), (80) 



With Eq. (80 ) in hand resolvent identities may be used to obtain the subtracted amplitude 



at any energy [20] : 

A<®(q, k; E) = A<®(q, k; -E D ) + B^(q, k; -E D ) + / dq'q' 2 Y^ \q, q'- E)A^{q', k; E) , 

Jo 

(81) 

where the second inhomogeneous term is given by 

B<%, k; E) = S[Z](q, k; E)+ t dq" q" 2 A<-°\q, q"; -E D ) r^(-E D ; q")5[Z](q", k; E) , (82) 

Jo 

with 

5[Z](q,q';E) = Z(q,q';E) - Z{q,q';-E D ) , (83) 



and t(E; q) = t(E - ^r). The kernel in Eq.(81) is 



Y(°\q, q'; E) = A< \q, q'; -E D )S[r^](E; q') + 5[Z](q, q'; E)r^(E; q') 

+ f A dq" q" 2 A^(q, q"; -E D )r^(-E D] q")S[Z](q" , q'; E)r^(E; q') 
Jo 

= A^(q,q';-E D )5[r^}(E;q') + B^(q,q';E)T^(E;q') , (84) 

with: 

5[r^mq') = r^(E;q')-r^(0;q'). (85) 

The result is a formulation of the three-body problem with short-range forces in which 
only renormalized quantities appear. Subtracted equations relate the amplitude at different 
energies (including energies below the atom-dimer scattering threshold). Hence, once 03 is 
known, the three-body bound-state spectrum and atom-dimer phase shifts can be computed. 

For pure S-wave interactions the only impact of effective-range corrections on the above 
argument is a modification of the two-body propagator r pU]. For arbitrary r s this can be 
written in the form 

2h 2 1 

T( ' ~ ~nM 2 -7 - ifc + f(7 2 + k 2 ) ' ( ' 

This two-body propagator has one additional pole which lies outside the region of validity 
of the EFT. Therefore, the two-body propagator has to be expanded up to the desired order 
in r s /a 4 : 

^ = iS' (87) 



Technically, r s must be replaced by p if we wish to consider any order beyond n = 1. 
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with for n < 3 

_ . 4 n 
i=0 

This form for is then inserted into the above equations and the calculation at NLO 
and N2LO can proceed exactly as with the LO calculation. Below we denote the resulting 
amplitude by A^ n \ 

B. Renormalization to 

When the subtracted integral equations are presented as has been done in Subsec- 



tion HI A it appears that to use them we must choose a 3 as the observable that provides the 
three-body input to the scattering equation. This constitutes a renormalization scheme that 
is appropriate for the low-energy properties of a system with a fixed atom-atom scattering 
length. In the situations of interest to us here though, the two-body scattering length of 
the atomic species is controlled by an external magnetic field, through the phenomenon of 
Feshbach resonances. The long-range properties of the atomic system (including the three- 
body scattering length) can then vary rapidly if the resonance is narrow. But, physically, 
the short-distance properties of the three-body system — in particular the asymptotic form of 
the three-body wave function in the high-momentum/short-distance region — should be un- 
affected by changes in the infra-red physics associated with the two-body scattering length. 
In what follows we will therefore use a renormalization scheme in which the parameter as- 
sociated with this ultra-violet part of the wave function is independent of the two-body 
scattering length. The analysis of Sec. [XT] shows that (in the unitary limit) the phase of this 
part of the wave function is in one-to-one correspondence with the binding momentum of the 
n*th three-body bound state, and so this renormalization scheme amounts to ensuring 
that the value of taken at the point 7 = 0, is unchanged when we include corrections 
that are proportional to r s . 

This means we need to know how to relate observables calculated in the unitary limit, 
7 = 0, to observables calculated at a finite value of 7. More generally, we would like 
to develop renormalized equations that encode relations between the amplitude A^ at 
different values of the two-body binding momentum 7. The 7-dependence of A^ arises 
in the integral equation through the dependence of r and Z on the two-body scattering 
length 5 . In what follows we use this fact, together with the subtraction approach, to relate 
the atom-dimer scattering amplitudes at two different values of 7. Choosing 7 = as one 
point in the resulting equation then allows us to compute A^ for any desired value of k*. 

Consider the operator form of the integral equation for the threshold scattering-amplitude 

A 

A(j, 0) = Z(y, 0) + Z(j, 0) r( 7 , 0) A(j, 0) , (89) 

where now we have written the 7-dependence explicitly, as the first argument in all quanti- 
ties, and indicated that this is an equation at threshold by including the atom-dimer relative 
momentum, k as a second argument. This supersedes the k = y/ME/h employed above, 



Z is 7-dependent at fixed on-shell momentum k, which is how we formulate our scattering problem. 
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and the argument E used there is now related to k by 



2 U2 



E 



3h z k 
AM 



E 



D- 



(90) 



At the atom-dimer scattering threshold k = 0. We have also suppressed the superscript n 
in Eq. (89) since the argument we develop in the next few paragraphs is independent of the 
T^ n > that is employed there. 

Equation (89) can be written as 

2(7 1 0) = i(7,0)[l + r( 7l 0)i( 7l 0)]~ 1 . (91) 



Assuming that we have already calculated the scattering amplitude for a given combination 
of 03 and a certain 7, say 71, using Eq. (77) and (80), we can now use Eq.(91) to derive the 
scattering amplitude for the two-body scattering length that corresponds to 72. Evaluating 
the expression for ^(71, 0) — 2 (72, 0) gives 



Z( 7l ,0)-Z( 72 ,0) = [1 + 4( 7 i,0)t( 7 



1, 1 



1-1 



.4(71,0) 



-*4( 72) 0)[1 + t( 72 , 0)4( T 2, 1 



1-1 



(92) 



Now multiplying with the appropriate expression from the left- and right-hand side, we 
obtain an integral equation for .4.(72, 0) 

4(72,0) = 4(7i, 0)+ 4( 7 i, 0)[r( 72 ,0)-r(7i,0)]4(72,0) 

+ [1 + 4( 7 i, 0)r( 7l , 0)] [Z( l2 , 0) - Z{ lu 0)] [1 + r( 72 , 0)4( 7 2, 0)] . (93) 

Once the threshold amplitude 4.(72, 0) is obtained from this equation, it can in turn be used 
to calculate the scattering amplitude at any energy and for the new value of 7, 72 . 
The strategy is thus: 

1. Pick a value of CZ3 and set 71 = 70 (chosen arbitrarily) and compute the threshold 4, 
4(7o,0), using Eqs. d77b and @. 



Use Eq. (93) to compute 4(0,0): the threshold amplitude in the unitary limit. 



Use Eq. (81) to obtain 4(7 = 0) for E < 0, thereby determining the bound-state 
spectrum. 



4. Adjust 03(70) and repeat steps 1-3 until a desired is obtained. 



5. Then use Eq. (93) to compute 4(7,0) for other values of 7. 



6. And once again use Eq. (81) to determine the bound-state spectrum. 

The results of this procedure are shown in Fig. [TJ This represents the spectrum of three- 
body bound states which correspond to different two-body scattering lengths but the same 
asymptotic behavior in the short-distance region (i.e. the same at 7 = 0). The binding 
energies are renormalized so that there is a state in the unitary limit which has binding 
momentum = VO.I570. (From now on we denote this state as the one with index n — 1.) 
We also set M — U — 1 . In Fig. [I] the binding momentum k — \/ MB3 /h is plotted against 
the inverse of the two-body scattering length. The solid line represents the leading-order 
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-3 -k 



FIG. 1: The bound-state spectrum of the three-boson system with short-range interactions. The 
solid lines denote the leading-order results and the dashed lines give the next-to-leading order 
result for an effective range of jor s = 0.01. Here the scale has been altered as per Ref. [I], with 
7 and n first rewritten in polar form as 7 = 70-ff cos£ and k = 70-ffsin^, and then the quantities 
H 1 / 4 cos £ and —if 1 / 4 sin ( displayed here. The dotted line denotes the point at which trimers 
become unstable against breakup to an atom and a dimer. 



results, which reproduce those presented in Ref. jTj. Although the axes are labelled — k 
and 7, we have actually plotted — if 1 / 4 sin £ versus if 1//4 cos£, with (the dimensionless) H 
and ( denned as 7 = 70-ffcosC and k = 70-ffsin^. This rescaling of the variables allows 
us to display a greater range 7 and k as the discrete scaling factor is reduced from 22.7 to 
22.7 1 / 4 = 2.2 pQ. 

We can extend these results to NLO by following Steps 1-7 above, but using the form (88) 
for n = 1 in the integral equations. We reiterate that we renormalize at NLO by demanding 
that the value of is unchanged from its LO value. Consequently the value of a 3 obtained 
at the starting point 7 = 70 undergoes a shift. In Fig. [I] we show an extended Efimov plot 
which includes the NLO results for a representative value of r s 7o = 0.01 as the dashed line. 
The effective range has been chosen such that the deepest three-body bound state shown 
remains in the range of validity of the EFT for all allowed values of 7. We see that the shift 
caused by range corrections is only recognizable for this state, and that even the shift in its 
binding energy vanishes in the unitary limit. 



C. Universality in the range correction 



We now take the results and methodology of the previous subsection and use them to 
compute corrections to the bound-state spectrum that are linear in the effective range r s . 
This calculation was performed in the limit \a\ — > 00 in Sec 
the case of finite a. In the small r s regime we can write (c.f. 



II B 



Eq. 



and here we extend it to 
21) 



2,^2 r 



B r . 



70 



7o 



(94) 
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FIG. 2: Left panel: The extracted function G\ for negative scattering length. Right panel: Results 
for G\ and rescaled G2 for positive scattering length The solid line denotes the result for the 
function G extracted from the bound state which has binding momentum = \/0.157o in the 
unitary limit. The dashed line denotes the rescaled result obtained from the bound state which 
has binding momentum k = k*/22.7 in the unitary limit. 



where the function F n gives the leading-order binding energy of the nth state at arbitrary 
7 pp. The scale 70 is arbitrarily chosen, and for the calculations presented here we have = 
\/0.157q. In general we will consider values of 7 which range from (i.e. the unitary limit) 
to 7 ~ 14. 2re*, where the n = 1 trimer state vanishes through the atom-dimer threshold. It 
needs to be noted, though, that at any non-zero 7 there are only a finite number of bound 
states within the domain of validity of the EFT: the Efimov spectrum of infinitely many 
shallow three-body bound states is only realized at 7 = 0. Indeed, at 7 = 70 there are only 
two three-body bound states which are inside the domain of validity of the EFT. 

Discrete scale invariance leads us to hypothesise that the G„'s are related to G n * by: 

where the factor out the front comes from assuming scaling of the first-order correction. 
If we computed solely in the unitary limit this pre-factor would be e ( n *- n ) 37r / s o ; a number 
of effects that modify this scaling for finite 7 (when the hyperspherical potential is no longer 
a sum of power laws, even for small r s ) are accounted for by the function 9 n . This obeys 

9 m {z) = 1 for all z, (96) 
n (O) = 1 for all n, (97) 



and is ~ 1 everywhere else. The constraint (97) on 9 n incorporates the result already 
demonstrated analytically in the limit 7 — > 0: that if the linear-range correction vanishes 
for one Efimov bound state then it vanishes for all of them. This result is verified by both 
our results and an independent computation [T8] . 



Here we choose n* — 1. The extraction of the function G n from Eq. (94) is numerically 
quite delicate. It involves a numerical derivative with respect to r s , but in order to keep 
consistency with the scales of our EFT we must maintain both r s 7 <C 1 and Ar s 3> 1. In 
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practice we compute G for a number of values of r s that satisfy these constraints and then 
extrapolate to the r s — > limit. 

The results for G\ for positive 7 are shown as the solid line in the right panel of Fig. [2j 
Away from 7 = the case of G does not give us a wide enough range of 7 over which 
to test Eq. (95) while still remaining inside the domain of validity of the EFT. We can, 
however, test the prediction (95) for the case n — 2. Since the binding energies for this 
case are two orders of magnitude smaller, significantly higher accuracy is needed to compute 
the necessary numerical derivatives. Results for the function G2 in the case of positive 
two-body scattering length are represented by the dashed line in the right panel of Fig. |2j 
with G2 having been rescaled according to Eq. (95). We see that the horizontal scale in 
the two cases is indeed related by a factor of e n < 1 °, exactly as predicted by discrete scale 
invariance. Also, both functions vanish at the same point, 7 ~ 5. 570, once this rescaling is 
performed (7 ~ 0.247 for n = 2 before rescaling). Moreover, the peak in the rescaled G2 
is at approximately the same position (7 pa 3. 570) as that in G±. But, although the factor 
of 22. 7 3 postulated in Eq. (95) allows functions that differ by four orders of magnitude to 
be plotted on the same vertical scale, the peak of G2 after rescaling is 30% higher than 
anticipated based on G\ and the scaling of binding energies in the unitary limit. Further 
investigation of these deviations from pure scaling that are encoded in the functions 8 n 
will be interesting, but are beyond the scope of this paper. 

That the rescaled G2 and the function Gi vanish at the same value of 7 would seem to be 
a simple consequence of the fact that the binding momentum of the three-body state goes 
to zero at this point, and so the first-order correction to the binding energy must be zero. 
But this straightforward observation has the important consequence that corrections linear 
in r s do not affect the values of 7 > at which states that would be part of the Efimov 
tower in the 7 = limit disappear from the bound-state spectrum. We conjecture that this 
result also applies to the case of negative scattering length, i.e. 7 < 0. 

In the left panel of Fig. [2] we display the result for G\ at negative scattering length. 
The rapid decrease of the three-body binding energies in this domain makes the extraction 
of the functions still more complicated, and in this case we achieve satisfactory numerical 
accuracy only for G\. The state with n — 1 vanishes above the three-body threshold at 
7 = — 0.26(l)7o, so we obtain G\ over almost the entire range where it is defined. The 
particularly fast variation of three-body binding energies near this vanishing point means 
that we could not extract the value of G\ at the endpoint with any serious precision. Thus 
our conjecture of the previous paragraph remains unconfirmed. But, if true, that conjecture 
has important implications for three-body recombination at negative scattering length. We 
therefore now turn our attention to how recombination can be computed from the scattering 
amplitudes we have already obtained. 



IV. RANGE CORRECTIONS TO THREE-BODY RECOMBINATION 

In this section we will consider effective-range corrections to the three-body recombination 
rate for positive two-body scattering length. The calculation of the rate at negative scat- 
tering length requires the evaluation of the S-matrix element for three-atom-to-three-atom 
scattering. This is more involved and will be considered in a later publication. 

Three-body recombination is a process in which three atoms collide to form a diatomic 
bound state and the energy that is thereby released causes the particles to leave the trap. 
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FIG. 3: The recombination rate a as a function of 7/70- The solid line gives the leading order 
result, while the dashed and dot-dashed lines give the results for r s 7o = 0.01 and 0.005, respectively. 

The rate of decrease in the number density of atoms in the trap is given by 

A n = -s an 3 . (98) 

The factor of 3 in the above equation arises from the assumption that all three particles 
involved in the recombination process will leave the trap. Atoms with large positive scat- 
tering length can recombine into the shallow dimer with Ejj = h 2/ y 2 /M. Recombination 
can also produce deep dimers with -Edcep ~ Ti 2 /(Ml 2 ) — if such states are supported by the 
underlying two-body interaction. For negative scattering length only recombination into 
such deep dimers is possible. 

Range corrections to three-body recombination have previously been considered in an 
EFT framework [21] . But there only correlations between the atom-dimer scattering length 
and a were considered. However, the effect that range corrections have on a as the two-body 
scattering length is varied but the three-body parameter (here denoted by k#) remains fixed 
is particularly relevant to experiment. This case corresponds to the situation in which the 
two-atom scattering length is tuned using an external magnetic field. Sizeable sensitivity 
of 7 (or equivalently a) to the magnetic field is called a Feshbach resonance and occurs if 
a change in the magnetic-field strength leads to a threshold crossing of a two-body bound 
state. 

Here, for simplicity, we consider three-body recombination at threshold. This corresponds 
to the zero-temperature limit. We will also neglect three-body recombination into deep 
dimers. An analytic form for the leading order recombination rate into shallow dimers at 
threshold has been derived by Macek et al. [22] and independently by Petrov [23] 



1287r 2 (47r -3V3)sin 2 (s ln^) % 
V ; sinh 2 (7rs ) + cos 2 (s ln^) M7 4 ' v ' 

The recombination rate therefore scales as the product of Ti/(M^ A ) and a log-periodic func- 
tion. The three-body parameter that determines the phase of this log periodicity defines 



22 



the scattering length at which the recombination rate has a minimum. At leading order it 
is related to the three-body parameter by a* = (e 7rs °) n 0.32fi:~ 1 []]. 

The recombination rate can be obtained from the elastic atom-dimer scattering amplitude 
A by solving the STM equation and evaluating [2TJ 121] 

647T 2 M 1 . , /- %l9 , v 

"= V3 (l- T V)'(l-7r) W0|2T/A0)l '- < 10 °> 



At leading order this reproduces Eq. (99). We now use the subtracted integral equation 
derived in Sec. Ill to calculate the amplitude A in Eq. (100). In any experiment the effective 
range will vary with the external magnetic field which is used to control the scattering length. 
However, this variation of r s will be slow compared to that of 7, and so we assume that the 
effective range stays constant as 7 varies. 

Figure [3] shows the resulting dependence of the threshold recombination rate a on 7 in 
units of The results display the log-periodic behavior expected from Eq. (100). 

The solid line gives the three-body recombination rate in the zero-range limit — the LO EFT 
result — while the dashed and dot-dashed lines give the NLO results for r s 7 = 0.01 and 
0.005, respectively. All calculations are renormalized so that they correspond to the same 
binding momentum « # of the first excited Efimov state at \a\ = 00. The renormalization at 
\a\ = 00 makes all curves approach the LO result as 7 — > 0. Away from 7 = it can be seen 
that the finite effective range leads to a shift of the minima in the recombination rate and 
also influences the maximum value of a. 

Given the recent experiment of Ref. [I] the extension of this computation to the case 
7 < will be particularly interesting. In particular, the minima in the recombination rate 
are determined by the values of the scattering length at which one of the Efimov states 
is at threshold. Therefore we conjecture that these minima receive no correction of 0(r s ). 
This will be discussed in a future publication. More generally for both 7 > and 7 < 
the inclusion of range corrections into calculations for the recombination of 133 Cs atoms will 
allow for a better description of data and therefore a more precise determination of the 
relevant three-body parameters. However, in order to obtain truly precise results it will be 
necessary to account both for the effects of deep dimers and the variation of the effective 
range as a function of the magnetic field. 



V. SUMMARY AND OUTLOOK 

Near a Feshbach resonance the scattering length, a, of the two-body system is much 
larger than all other length scales in the problem. This leads to an Efimov spectrum in the 
three-body system, which displays discrete scale invariance. However, corrections due to 
the effective range, r s , will in principle always affect observables even in the limit of infinite 
scattering length. This paper considered the impact of such corrections on three-body system 
observables. 

We did this first by examining the spectrum of the three-body system in the unitary limit, 
where \a\ — > 00. Using hyperspherical coordinates for the three-body system, we extended 
previous analyses [HE5] to derive the potential that corrects the strict 1/R 2 potential present 
in the unitary limit for the presence of a non-zero r s . (This potential had been conjectured, 
but not derived in Ref. [IT].) Specifically, we found that the perturbing potential is 

W(R) = -^t^, (101) 
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where so ~ 1.00624 is the eigenvalue of the leading-order hyper-angular equation and £o ~ 
0.480 is obtained by analyzing that equation for r s ^ 0. 

We assessed the impact of the potential on the Efimov spectrum of three-body 
bound states. Since the potential is singular renormalization is required. If we take as our 
renormalization condition that the n*th state is unperturbed by the correction, then the 
shift in the nth state is zero in the limit A — * oo, and so the entire spectrum is unperturbed. 
This is consistent with the discrete scale invariance that is present in that limit, but we note 
that the conclusion is dependent on the particular renormalization condition chosen. 

We then examined how these results change as we move away from the limit \a\ — > oo. We 
derived subtracted integral equations that allowed us to relate the 1 + 2 scattering amplitude 
at different values of the scattering length a. This led us to an extended version of the Efimov 
plot, where we were able to display how three-body bound state energies vary not only with 
a, but also with r s . Such a plot could also be generated with an explicit three-body force, but 
the subtractive approach employed here has certain advantages, in particular when higher- 
order calculations at large cutoffs are pursued. This approach has been used with success in 
Refs. pm EED, EU [26] in a variety of nuclear and atomic-physics contexts to compute results 
that are claimed to be valid up to N2LO in the r s /a expansion. It is important to point 
out that the disagreement between these works and Refs. [9j [121 123 [28] regarding whether 
additional three-body input is needed to renormalize the three-body problem at N2LO does 
not affect the results presented here. The formalism employed to compute the first-order 
correction in r s reproduces the results of Refs. [HJ H21 123 EE] , even if it leads to conclusions 
at N2LO that disagree with these studies. 

Thus the only assumption in our work is that \r s \ ^> /, where I is the underlying length 
scale in the two-body potential. Because of this assumption all calculations of Sec. Ill 
were performed in the limit A — > oo.We extracted, in this limit, a universal function that 
describes the linear (in r s ) correction to the bound-state energies of the three-body bound- 
state spectrum. We chose to examine the correction for the state with n — n* — 1, with the 
universal function defined by: 



2,^2 r 



B 



" ~ ~M 



K* ( — )+«* r s G m ( — ) + 0[(«» r s 



(102) 



where F n is given in Ref . [T] , and defines the 7 dependence of the bound state in the r s = 
(scaling) limit. The result we obtained for G n * is displayed in Fig. [2j We have furthermore 
showed that the next-to-leading order corrections to different states in the Efimov spectrum 



are approximately related to each other by the scale transformation defined in Eq. (95). It 
is an interesting question whether this symmetry of the leading-order wave function impacts 
other observables too. 

Being able to use subtracted equations to relate physical observables at different scat- 
tering lengths to each other led us to consider three-body recombination into the shallow 
dimer. We calculated the recombination rate for different values of the effective range r s 
while renormalizing to the binding momentum of an Efimov state in the unitary limit. 
Consequently our results for different values of the two-body effective range approach one 
another as the two-body scattering length increases. 

In order to examine the behavior of the three-body system near a Feshbach resonance 
we have used an EFT that exploits the hierarchy of scales \r s \ <C \a\. We have focused on 
effects of 0(r s /a). This enables us to probe the pattern of convergence of the EFT, thereby 
facilitating reliable error estimates for observables calculated in this framework. It also 
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provides a clear differentiation between effects that are universal (i.e. arise solely as a result 
of \a\ ^> \r s \) and effects that depend on details of the underlying two-body interaction. It is 
interesting to speculate that features such as the Phillips and Tjon lines in nuclear physics 
are well reproduced at leading order [23 ED], in spite of the large expansion parameter there 
{r s /a ps 0.3) because the discrete scale invariance suppresses the NLO correction below the 
natural expectation. 

However, in order to confirm this speculation an analysis of the four-body system in, 
and near, the unitary limit must be carried out. Such a conjecture might also be rendered 
more plausible if the discrete scale invariance that is present at leading order suppressed 
some of the next-to-next-to-leading order corrections to the Efimov spectrum. In order to 
examine this possibility we would like to extend the above analysis to 0(r 2 /a 2 ). While 
the computation of such corrections to the bound-state spectrum appears straightforward, 
it is numerically somewhat delicate. However, a computation of the 0(r 2 ) corrections to 
the unitary-limit spectrum can be carried out using the hyperradial formalism in a manner 
similar to that employed here. 
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APPENDIX A: A PROBLEM WITH r s7 < 

For r s 7 < problems with the method we have adopted to obtain the correction ~ r s in 
the bound-state spectrum occur as soon as we consider cutoffs A 3> 1/r. In this case the 



largest eigenvalues of the kernel of the original STM equation (72) are negative, and so the 



subtraction does not improve the equation's behavior. Instead, when we solve the subtracted 



integral equation (81) at NLO we obtain spurious low-energy solutions which have no LO 
counterpart, i.e. they are not smoothly connected to a LO eigenvector by variation of r s . 
Examining these eigenvectors we see that they have most of their support at p ~ 1/r. These 
eigenvectors are thus a non-perturbative effect generated by the inclusion of the pieces ~ r s 
in the kernel. The effect of these NLO pieces of the kernel should be perturbatively small, 
so the fact that they produce eigenvectors with support at p ~ 1/r is not a valid prediction 
of the EFT. 

Therefore we wish to eliminate these eigenvectors from the spectrum of the kernel, leaving 
us with only eigenvectors (and eigenvalues) which are perturbatively close to those that exist 
at leading order. We can do this by noting that, regardless of the sign of r s , the integral 
equation with the N2LO kernel above yields small changes from the LO result for small r s . 



The largest eigenvalues of the N2LO kernel (n = 2 in Eq. (88)) are positive regardless of 
the sign of r s — as in the LO case. Thus the subtraction approach works straightforwardly 
there. Therefore, we choose to stabilize the troublesome NLO integral equation for NLO by 
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adding an admixture of the N2L0 kernel to the NLO one. In other words, we choose for 
the residue function S(E; q) 



2fr 4 / 3 , mE \ . .3 



Sri, o , os / /3 „ ME 



■ ft ¥ (p2 " fc2) 7 i/"i 11 ^ (A1) 



with i£ and related by E = — -Ed- Here the piece on the second line effects the 

stabilization, but note that it is not the full N2LO correction, since we have omitted pieces 
of the kernel that scale as (7r s ) 2 , as they do not assist with the stabiliz ation . Therefore 
Si S^ 2 \ This distinction is irrelevant in the unitary limit, and so there (Al) interpolates 



smoothly between the NLO result (for a = 0) and the N2LO result for (a = 1). And, 
provided perturbation theory applies, the result of using Sjp instead of should be linear 
in a. Therefore as long as r s is small enough that perturbation theory is valid we can 
examine the results for calculations with S^ 14 ^ and make a linear extrapolation to a = in 
order to obtain the NLO result. 

Below we show a table of results for B\ using for r s 7 = —0.07, and r s 7 = —0.14 
and 037 = 1.79. Once a < 0.3 the admixture of the N2LO piece in the propagator no longer 
stabilizes the integral equation. However, the data for a > 0.4 is already sufficient to allow 
us to predict: 

B < i ) lE D {r sl = -0.07) = 1.7115(1), (A2) 

and 

B l i ) /E D {r sl = -0.14) = 1.6983(1), (A3) 
for the pure NLO piece of the result. When combined with the LO result: 

B ( i ) /E D {r sl = 0.0) = 1.723, (A4) 

this provides nice evidence of linearity in r s 7, even though we are now considering r s < 0. 
Moreover, these predictions are entirely consistent with the NLO results for positive r^: 

Bi 1} /E D (r sl = 0.07) = 1.734; B { i ] /E D {r sl = 0.14) = 1.743, (A5) 

and the assumption of linearity in r^. 

Indeed, ultimately this is the solution to the instability we are attempting to cure: if 
perturbation theory is valid then the shift at NLO should be linear in r s . Since we can 
compute it for positive r s using the methods outlined above we can extrapolate those results 
linearly to the region r s < and thereby obtain the perturbative shift in bound-state 
energies. If that result is different from what is obtained by the use of the NLO kernel with 
r s < then that difference is an effect beyond perturbation theory in r s 7 and r s K, and such 
effects are not our concern here. 
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